In [1]:
%matplotlib inline
%config InlineBackend.figure_format='retina'
from IPython.display import Image, HTML
import random
import math
import pysolar
import datetime as dt
import pytz
import matplotlib.dates as mpd
import scipy
import numexpr as ne
from datacharm import *
from datacharm.dataplots import tableau20, ch_color
from enerpi.base import timeit
from enerpi.api import enerpi_data_catalog
from enerpi.enerplot import plot_tile_last_24h, plot_power_consumption_hourly
from enerpi.process import separa_ldr_artificial_natural, get_solar_day
from enerpi.radsolar import ajuste_ldr_con_barrido_irrad_inc, local_maxima_minima
LAT = 38.631463
LONG = -0.866402
ELEV_M = 500
TZ = pytz.timezone('Europe/Madrid')
FS = (16, 10)
def _tslim(ax, h0=12, hf=22, m0=0, mf=0):
xlim = [mpd.num2date(x, tz=TZ) for x in ax.get_xlim()]
new = [mpd.date2num(d.replace(hour=h, minute=m)) for d, h, m in zip(xlim, [h0, hf], [m0, mf])]
ax.set_xlim(new)
return ax
def set_sun_times(df, tz=TZ, lat_deg=LAT, long_deg=LONG, offset='5min'):
if df.index[-1].date() == df.index[0].date():
# Same day:
ts_day = pd.Timestamp(df.index[0].date(), tz=tz)
sunrise, sunset = pysolar.util.get_sunrise_sunset(lat_deg, long_deg, ts_day + pd.Timedelta('12h'))
delta = pd.Timedelta(offset)
sunrise, sunset = [sunrise - delta, sunset + delta]
df['sun'] = False
df.loc[sunrise:sunset, 'sun'] = True
return df
else:
print_err('IMPLEMENTAR!')
assert()
sns.set_style('whitegrid')
# Catálogo y lectura de todos los datos.
cat = enerpi_data_catalog()
data, data_s = cat.get_all_data()
print_info(data_s.tail())
LDR = pd.DataFrame(data.ldr).tz_localize(TZ)
print_cyan(LDR.describe().T.astype(int))
LDR.hist(bins=(LDR.max() - LDR.min() + 1).values[0] // 5, log=True, figsize=(18, 8))
Image('/Users/uge/Desktop/LDR_analysis/LDR_análisis_día_2016_09_12.png')
Out[1]:
In [2]:
#day = '2016-08-27'
day = '2016-09-12'
data_raw = LDR.loc[day:day]
print_info(data_raw.head(3))
data, data_simple = separa_ldr_artificial_natural(data_raw, resample_inicial='2s', delta_roll_threshold=100)
solar = get_solar_day(data)
cols_p = ['ldr', 'delta_roll','artif_level_max', 'altitud','irradiacion_cs']
solar[cols_p].iloc[::5].plot(figsize=FS)
Out[2]:
In [8]:
from scipy.signal import medfilt
#panel = pd.read_hdf('/Users/uge/Desktop/LDR_analysis/POIS.h5', 'POIS')
#panel.loc['2016-08-28']
#ts = pd.Timestamp('2016-08-28')
#df = panel.loc['2016-08-28']
#df['day'] = ts.date()
#from enerpi.radsolar import identifica_pois_ldr
@timeit('identifica_pois_LDR', verbose=True)
def _identifica_pois_ldr(data_raw):
def _es_max_local(x):
mid = len(x) // 2
if (np.argmax(x) == mid) and (x[mid] > 10):
return True
return False
def _pendiente(x):
return (x[-1] - x[0]) / len(x)
def _pendiente_ini(x):
if len(x) > 100:
return _pendiente(x.values[1:101])
return _pendiente(x.values)
def _pendiente_fin(x):
if len(x) > 100:
return _pendiente(x.values[-100:-1])
return _pendiente(x.values)
def _pendiente_last_gen(x):
limit = min(len(x), 600)
if limit > 1:
return _pendiente(x.values[-limit:-1])
print_err(x)
return np.nan
def _get_alt(d):
alt = pysolar.solar.get_altitude_fast(LAT, LONG, d)
if alt > 0:
return alt
return 0
# Resampling 1s + sun
data_homog = data_raw.resample('1s').mean().interpolate().dropna()
data_homog = set_sun_times(data_homog, tz=TZ, lat_deg=LAT, long_deg=LONG, offset='5min')
subset_pysolar = data_homog.iloc[::120].index
data_homog = data_homog.join(pd.Series(subset_pysolar, index=subset_pysolar, name='altitud'
).apply(_get_alt)).interpolate()
# Median filter
data_homog['median_filter'] = medfilt(data_homog.ldr, kernel_size=301)
# BUSCA POI's
data_homog['r_disp'] = data_homog.median_filter.rolling(5, center=True).apply(lambda x: x.max() - x.min()).fillna(0)
data_homog['poi'] = data_homog['r_disp'].rolling(5).apply(_es_max_local).fillna(False).astype(bool)
# Construye intervalos entre POI's
data_homog['intervalo'] = data_homog['poi'].cumsum()
# DATOS POIS:
gb_poi = data_homog.tz_localize(None).reset_index().groupby('intervalo')
pois = pd.DataFrame(pd.concat([gb_poi.sun.apply(lambda x: np.sum(x) / len(x)).rename('fr_sun'),
gb_poi.ldr.count().rename('seconds'),
gb_poi.altitud.mean(),
gb_poi.median_filter.apply(lambda x: (x.values[-1] - x.values[0]) / len(x)
).rename('pendiente_tot'),
gb_poi.median_filter.apply(lambda x: np.round(np.median(x.values[1:10]))
).rename('inicio'),
gb_poi.median_filter.apply(lambda x: np.round(np.median(x.values[-10:-1]))
).rename('fin'),
gb_poi.median_filter.apply(_pendiente_ini).rename('pendiente_i'),
gb_poi.median_filter.apply(_pendiente_fin).rename('pendiente_f'),
gb_poi.median_filter.apply(_pendiente_last_gen).rename('pendiente_last_mf'),
gb_poi.ts.first().rename('ts_ini'),
gb_poi.ts.last().rename('ts_fin')], axis=1))
pois['step'] = (pois['inicio'] - pois['fin'].shift()) # .fillna(0)
print_ok(pois.head(10))
return data_homog, pois
df = LDR.loc['2016-08-15']
data_homog, pois = _identifica_pois_ldr(df)
pois
Out[8]:
In [10]:
data_homog.ldr.plot()
Out[10]:
In [26]:
def _plot_bt_medfilt(data_homog, t0='0:00', tf='1:00', ks_1=301, ks_2=75):
df = data_homog.ldr.between_time(t0, tf)
ax = df.plot(figsize=FS, color=tableau20[6], lw=.75, alpha=.7)
pd.Series(medfilt(df, kernel_size=ks_1), index=df.index,
name='Medfilt_301').plot(ax=ax, color=tableau20[4], lw=4, alpha=.5)
pd.Series(medfilt(df, kernel_size=ks_2), index=df.index,
name='Medfilt_301').plot(ax=ax, color=tableau20[2], lw=3, alpha=.6)
return ax
_plot_bt_medfilt(data_homog, t0='0:00', tf='1:00', ks_1=301, ks_2=75)
plt.show()
_plot_bt_medfilt(data_homog, t0='15:00', tf='18:00', ks_1=301, ks_2=75)
plt.show()
_plot_bt_medfilt(data_homog, t0='19:00', tf='23:59:59', ks_1=301, ks_2=75);
In [76]:
@timeit('SUNTIMES', verbose=True)
def _suntimes(df):
tt = pd.DatetimeIndex(df.index.date, tz=TZ)
delta = pd.Timedelta('10min')
s_sun = pd.Series([False] * len(df), index=df.index, name='sun')
for ts_day in pd.DatetimeIndex(set(df.index.date), tz=TZ):
sunrise, sunset = pysolar.util.get_sunrise_sunset(LAT, LONG, ts_day + pd.Timedelta('12h'))
sunrise, sunset = [sunrise - delta, sunset + delta]
s_sun.loc[sunrise:sunset] = True
data_ss = LDR.resample('5s').max()
df = _suntimes(data_ss)
In [203]:
#pois = pd.read_hdf('/Users/uge/Desktop/LDR_analysis/POIS_debug.h5', 'POIS').iloc[:-1].fillna(0)
#print_info(pois.describe())
sns.pairplot(pois_labeled.drop(['ts_ini', 'is_poi_to_label'], axis=1))
#pois.drop(['ts_ini', 'ts_fin'], axis=1)
Out[203]:
In [250]:
from datacharm.datafuncs import trimedia_TM, kmeans_1d, kmeans_md, clustering_pca_2d, silhouette_analysis, KMeans, PCA
X = pois_labeled.drop(['ts_ini', 'is_poi_to_label', 'intervalo'], axis=1).astype(float).values
#centers_km_c, labels_km_c, X_km_c, km_c, df_cl_c = kmeans_1d(pois_labeled, 'step', 7) #, func_preprocess=np.sqrt)
#silhouette_analysis(centers_km_c, labels_km_c, X_km_c)
centers_km_c, labels_km_c, X_km_c, km_c = kmeans_md(pois_labeled, ['step', 'altitud', 'ldr_median', 'pendiente_tot'], 4, sort_labels=False)
silhouette_analysis(centers_km_c, labels_km_c, X_km_c)
In [247]:
# pois.columns
'''['intervalo', 'fr_sun', 'ini_sun', 'fin_sun', 'seconds', 'altitud', 'ldr_mean', 'ldr_median', 'pendiente_tot',
'inicio', 'fin', 'pendiente_i', 'pendiente_f', 'pendiente_last_mf', 'ts_ini', 'step', 'is_poi_to_label',
'solo_natural', 'artificial_princ']'''
cols_clustering = ['fr_sun', 'ini_sun', 'fin_sun', 'seconds', 'altitud', 'ldr_mean', 'ldr_median', 'pendiente_tot',
'inicio', 'fin', 'pendiente_i', 'pendiente_f', 'pendiente_last_mf', 'step']
centers, labels, X, km = kmeans_md(pois_labeled, cols_clustering, k, sort_labels=False)
In [248]:
len(pois), len(pois_labeled), len(pois_to_calc)
Out[248]:
In [251]:
pois_km = pois.drop(['ts_ini', 'intervalo'], axis=1).copy()
#pois_km['km_group'] = labels
pois_km.head()
Out[251]:
In [291]:
#pois_km[pois_km.artificial_princ.notnull() & pois_km.artificial_princ]
#pois.drop(['ts_ini', 'artificial_princ', 'solo_natural', 'is_poi_to_label', 'intervalo'], axis=1)
In [293]:
#saltos = pois.join((pois['fin'] - pois['inicio'].shift(-1)).rename('step_salto'))
#saltos = saltos.join((saltos['pendiente_f'] - saltos['pendiente_i'].shift(-1)).rename('pendiente_salto'))
#saltos = saltos.iloc[1:-1][['seconds', 'fin_sun', 'fin', 'solo_natural', 'artificial_princ', 'step_salto', 'pendiente_salto']]
#saltos_labeled = saltos[saltos.solo_natural.notnull()]
#saltos_to_calc = saltos[saltos.solo_natural.isnull()]
#ts = saltos_to_calc.index.time[9]
#saltos.seconds.hist(bins=1000, figsize=FS, lw=0, log=True)
#(saltos.seconds).value_counts().sort_index()
In [262]:
k = 9
km_pca = clustering_pca_2d(X, k, whiten=True, ax=None)
X_2 = pois.drop(['ts_ini', 'artificial_princ', 'solo_natural', 'is_poi_to_label', 'intervalo'], axis=1).astype(float).iloc[1:].values
km_pca_2 = clustering_pca_2d(X_2, k, whiten=True, ax=None)
In [316]:
signal.general_gaussian?
In [346]:
from scipy import signal
#pd.set_option('display.width', 120)
path_st = '/Users/uge/Desktop/LDR_analysis/POIS_ident_1s.h5'
# data_homog, data_rs, pois = genera_pois_de_all_data_LDR(path_st)
pois = pd.read_hdf(path_st, 'POIS')
data_homog = pd.read_hdf(path_st, 'data_homog')
data_rs = pd.read_hdf(path_st, 'data_rs')
#pois_labeled = pois[pois.solo_natural.notnull()]
#pois_to_calc = pois[pois.solo_natural.isnull()]
#print(len(pois), len(pois_labeled), len(pois_to_calc))
print(len(pois), pois.columns)
data = data_rs.join(data_homog.ldr.rename('raw_ldr')).drop('sun', axis=1)
#data.reset_index().iloc[50:100]
window = signal.general_gaussian(51, p=0.5, sig=20)
filtered = signal.fftconvolve(window, data.raw_ldr)
filtered = (np.average(data.raw_ldr) / np.average(filtered)) * filtered
filtered = np.roll(filtered, -25)
print(len(filtered), len(data), window)
data['filter_g'] = filtered[:-50]
data.iloc[6000:8000].drop(['intervalo', 'altitud'], axis=1).plot(figsize=FS, lw=2, alpha=.7)
Out[346]:
In [375]:
import statsmodels.api as sm
#sm.tsa.ARIMA?
day = data.loc['2016-09-04']
print(day.shape)
#data.median_filter.shape
day.plot()
#data.median_filter.rolling(10).apply(np.ptp)
#np.ptp?
Out[375]:
In [362]:
data['mf_dmax'] = data['median_filter'].rolling(10).max() - data['median_filter']
data['mf_dmin'] = data['median_filter'] - data['median_filter'].rolling(10).min()
data['mf_dif'] = data['mf_dmax'] - data['mf_dmin']
#data['median_filter'].rolling(10).max() - data['median_filter'].rolling(10).min()
ax = data[['mf_dif', 'median_filter']].iloc[5000:7000].plot(figsize=FS, lw=2, alpha=.6) #, color=tableau20[0])
#data['median_filter'].rolling(10).min().iloc[6000:7000].plot(ax=ax, color=tableau20[2])
#data['median_filter'].rolling(10).max().iloc[6000:7000].plot(ax=ax, color=tableau20[4])
#
In [365]:
#data[data['mf_dif'].abs() > 15]
data['mf_dif']
Out[365]:
In [304]:
saltos = pois[['seconds', 'alt_sun_fin', 'fin', 'pendiente_f', 'pendiente_last_mf',
'ts_fin', 'step', 'pendiente_salto', 'is_poi_to_label']].copy()
saltos
Out[304]:
In [172]:
ax = ldr_raw.iloc[::5].plot(figsize=FS, lw=1, color='darkorange', alpha=.9)
(dplot['altitud'] * 12).plot(ax=ax, lw=4, color='y', alpha=.6)
dplot['median_filter'].plot(ax=ax, color='darkred', lw=2, alpha=.7)
Out[172]:
In [178]:
from pykeyboard import PyKeyboard, mac
k = PyKeyboard()
k.press_key('command')
k.tap_key('tab')
k.release_key('command')
In [54]:
#data_homog.ldr.between_time().rolling(75, center=True, win_type='triang').sum().plot()
t0, tf = '0:00', '1:00'
ks_1, ks_2 = 301, 75
df = data_homog.ldr.between_time(t0, tf)
#ax = df.plot(figsize=FS, color=tableau20[6], lw=.75, alpha=.7)
s1 = pd.Series(medfilt(df, kernel_size=ks_1), index=df.index,
name='Medfilt_301').resample('5s').median() #.plot(ax=ax, color=tableau20[4], lw=4, alpha=.5)
s2 = pd.Series(medfilt(df, kernel_size=ks_2), index=df.index,
name='Medfilt_301').resample('5s').median() #.plot(ax=ax, color=tableau20[2], lw=3, alpha=.6)
s3 = pd.Series(medfilt(df, kernel_size=11), index=df.index,
name='Medfilt_301').resample('5s').median() #.plot(ax=ax, color=tableau20[8], lw=2, alpha=.6)
ax = (s2 - s1).plot(figsize=FS, color=tableau20[6], lw=.75, alpha=.7)
(s3 - s2).plot(ax=ax, color=tableau20[4], lw=1.5, alpha=.7)
Out[54]:
In [3]:
def _es_max_local(x):
mid = len(x) // 2
if (np.argmax(x) == mid) and (x[mid] > 10):
return True
return False
def _pendiente(x):
return (x[-1] - x[0]) / len(x)
def _pendiente_ini(x):
if len(x) > 100:
return _pendiente(x.values[1:101])
return _pendiente(x.values)
def _pendiente_fin(x):
if len(x) > 100:
return _pendiente(x.values[-100:-1])
return _pendiente(x.values)
def _pendiente_last_gen(x):
limit = min(len(x), 600)
return _pendiente(x.values[-limit:-1])
def _get_alt(d):
alt = pysolar.solar.get_altitude_fast(LAT, LONG, d)
if alt > 0:
return alt
return 0
# Resampling 1s + sun
data_homog = data_raw.resample('1s').mean().interpolate().dropna()
data_homog = set_sun_times(data_homog, tz=TZ, lat_deg=LAT, long_deg=LONG, offset='5min')
subset_pysolar = data_homog.iloc[::120].index
data_homog = data_homog.join(pd.Series(subset_pysolar, index=subset_pysolar, name='altitud'
).apply(_get_alt)).interpolate()
# Median filter
data_homog['median_filter'] = scipy.signal.medfilt(data_homog.ldr, 301)
# BUSCA POI's
data_homog['r_disp'] = data_homog.median_filter.rolling(5, center=True).apply(lambda x: x.max() - x.min()).fillna(0)
data_homog['poi'] = data_homog['r_disp'].rolling(5).apply(_es_max_local).fillna(False).astype(bool)
# Construye intervalos entre POI's
data_homog['intervalo'] = data_homog['poi'].cumsum()
# DATOS POIS:
gb_poi = data_homog.tz_localize(None).reset_index().groupby('intervalo')
pois = pd.concat([gb_poi.sun.apply(lambda x: np.sum(x) / len(x)).rename('fr_sun'),
gb_poi.ldr.count().rename('seconds'),
gb_poi.altitud.mean(),
gb_poi.median_filter.apply(lambda x: (x.values[-1] - x.values[0]) / len(x)).rename('pendiente_tot'),
gb_poi.median_filter.apply(lambda x: np.round(np.median(x.values[1:10]))).rename('inicio'),
gb_poi.median_filter.apply(lambda x: np.round(np.median(x.values[-10:-1]))).rename('fin'),
gb_poi.median_filter.apply(_pendiente_ini).rename('pendiente_i'),
gb_poi.median_filter.apply(_pendiente_fin).rename('pendiente_f'),
gb_poi.median_filter.apply(_pendiente_last_gen).rename('pendiente_last_mf'),
gb_poi.ts.first().rename('ts_ini'),
gb_poi.ts.last().rename('ts_fin')], axis=1)
pois['step'] = (pois['inicio'] - pois['fin'].shift()) # .fillna(0)
pois.head(10)
Out[3]:
In [4]:
# PLOT
def plot_data_ldr(data_homog):
dplot = data_homog.tz_localize(None)
ax = dplot.ldr.iloc[::5].plot(figsize=FS, lw=1, color='darkorange', alpha=.9)
(dplot['altitud'].iloc[::10] * 12).plot(ax=ax, lw=4, color='y', alpha=.6)
dplot['median_filter'].iloc[::5].plot(ax=ax, color='darkred', lw=2, alpha=.7)
(dplot['intervalo'] * 20).plot(ax=ax, lw=3, color='magenta', alpha=.6)
ylim = ax.get_ylim()[1]
ax.vlines(dplot[dplot['poi']].index, 0, ylim, lw=2, color='darkgrey', alpha=.8, linestyle='--')
ax.xaxis.set_major_locator(mpd.HourLocator())
ax.xaxis.set_major_formatter(mpd.DateFormatter('%H:%M'))
return ax
ax = plot_data_ldr(data_homog)
ylim = ax.get_ylim()[1]
intervalos_largos = pois[pois.seconds > 600].copy()
intervalos_largos['ts'] = intervalos_largos['ts_ini'] + (intervalos_largos['ts_fin'] - intervalos_largos['ts_ini']) / 2
for i, t in intervalos_largos['ts'].items():
ax.annotate(str(i), (t, ylim - 50), color='red', ha='center')
plt.show()
In [ ]:
In [7]:
#momentos_sol = gb_poi.sun.apply(lambda x: np.sum(x)).rename('fr_sun').cumsum()
#(momentos_sol.max() - momentos_sol) / momentos_sol.max()
intervalo_luz_natural = []
intervalo_luz_artificial_principal = []
offset_step = pd.Timedelta('5min')
d_calc = data_homog.tz_localize(None)
ymax = ((d_calc.median_filter.max() + 100) // 100) * 100
for interv, row in pois.iterrows():
print_secc('Intervalo {}. De {:%d/%m/%y %H:%M:%S} a {:%H:%M:%S} ({} segundos)'
.format(interv, row.ts_ini, row.ts_fin, row.seconds))
delta_interv = row.ts_fin - row.ts_ini
offset = delta_interv / 10
df_slice = d_calc.loc[row.ts_ini - offset: row.ts_fin + offset]
ax = plot_data_ldr(df_slice)
#ax =plt.gca()
ini_p = mpd.date2num(row.ts_ini - offset)
fin_p = mpd.date2num(row.ts_fin + offset)
print(fin_p - ini_p)
#ax.add_patch(mpp.Rectangle((ini_p, 0), fin_p - ini_p, ymax, angle=0.0, color='red', linewidth=100, fill=True))
#ax.fill_betweenx([0, ymax], mpd.date2num(row.ts_ini - offset), mpd.date2num(row.ts_fin + offset), alpha=.6, color=tableau20[1])
#ax.fill_betweenx([0, ymax], row.ts_ini - offset, row.ts_fin + offset, alpha=.6, color=tableau20[1])
ax.fill_betweenx([0, ymax], row.ts_fin - offset_step, row.ts_fin + offset_step, alpha=.4, color=tableau20[1])
#ax.autoscale()
break
#dplot
In [9]:
def _str_to_bool(resp):
resp = resp.lower().rstrip().lstrip()
if resp in ['si', 'sí', 's', '1', 'yes', 'y']:
return True
return False
def _get_input():
es_nat = _str_to_bool(input('Es intervalo natural? '))
es_art = _str_to_bool(input('Es cambio por luz artificial? '))
return es_nat, es_art
_get_input()
Out[9]:
In [290]:
pois[pois.pendiente_last_mf.abs() < 1]
Out[290]:
In [196]:
def _rs_raw(x):
n = len(x)
return np.round(np.sum(x) / n**2, 0) * n
#data_homog.rolling(5, center=True, min_periods=1).apply(_rs_raw).plot(figsize=FS)
In [197]:
# pysolar.util.diffuse_underclear(latitude_deg, longitude_deg, when, elevation=0.0, temperature=288.15,
In [7]:
def _calc_pendiente(serie_homogenea, roll_window=3):
return serie_homogenea.diff().rolling(roll_window, center=True
).apply(lambda x: (x[-1] - x[0]) / roll_window
).rename('{}_roll_{}'.format(serie_homogenea.name, roll_window))
pd.concat([_calc_pendiente(data_homog.ldr, 5),
_calc_pendiente(data_homog.ldr, 30),
_calc_pendiente(data_homog.ldr, 120)], axis=1).plot(figsize=FS, alpha=.5, lw=1)
Out[7]:
In [8]:
ts_day = pd.Timestamp(data_homog.index[0].date(), tz=TZ)
print_yellowb(ts_day)
#pysolar.util.direct_underclear?
# NO -> pysolar.radiation.get_radiation_direct?
data_calc = pd.concat([data_homog,
_calc_pendiente(data_homog.ldr, 5),
_calc_pendiente(data_homog.ldr, 120),
_calc_pendiente(data_homog.ldr, 6000)], axis=1)
#.plot(figsize=FS, alpha=.5, lw=1)
#data_calc['dif_30'] = data_calc['ldr_roll_5'] - data_calc['ldr_roll_30']
data_calc['dif_120'] = data_calc['ldr_roll_5'] - data_calc['ldr_roll_120']
data_calc['dif_6000'] = data_calc['ldr_roll_5'] - data_calc['ldr_roll_6000']
#ax = data_homog[~data_homog.sun].plot(figsize=FS, lw=1, color='violet')
#data_homog[data_homog.sun].plot(ax=ax, lw=2, color='darkorange')
#data_homog[data_homog.sun].plot(figsize=FS, lw=2, color='darkorange')
data_calc[data_calc.sun]['ldr_roll_6000'].plot(figsize=FS)
print(data_calc.head())
data_calc.iloc[:-25000:1000].tail(10)
Out[8]:
In [ ]:
In [ ]:
In [9]:
rs = data_raw.ldr.resample('3s')
data_3 = pd.concat([rs.max().rename('ldr_max'),
rs.min().rename('ldr_min'),
rs.mean().rename('ldr_mean'),
rs.median().rename('ldr_median')], axis=1)
data_3['pend'] = rs.apply(lambda x: (x[-1] - x[0]))
data_3['disp'] = data_3['ldr_max'] - data_3['ldr_min']
data_3 = set_sun_times(data_3)
data_3.head(20)
Out[9]:
In [198]:
#data_3['rms'] = rs.apply(lambda x: np.sqrt(np.mean(x**2)))
#data_3['rms_roll'] = data_3.rms.rolling(10, center=True).apply(lambda x: np.sqrt(np.mean(x**2)))
#data_3['delta_rms'] = data_3['rms'] - data_3['rms_roll']
In [ ]:
In [ ]:
In [36]:
delta = pd.Timedelta('10s')
for t in maximos.index:
delta_ml = maximos.loc[t]
y_ldr = data_homog.loc[t, 'ldr']
y_mf = data_homog.loc[t, 'median_filter']
with_sun = data_homog.loc[t, 'sun']
mf_slice = data_homog.loc[t - delta:t + delta, ['ldr', 'median_filter']]
print_info('{} -> {} (delta_ml={}, median_filter={}, sun={})\n{}'.format(t, y_ldr, delta_ml, y_mf, with_sun, mf_slice))
#plt.plot(t, , marker='o', color='red', markersize=5)
break
maximos
Out[36]:
In [73]:
ax = data_homog.ldr.iloc[::5].plot(figsize=FS, lw=1, color='darkorange', alpha=.9)
data_homog['median_filter'].iloc[::5].plot(ax=ax, color='darkred', lw=2, alpha=.7)
for t in maximos.index:
plt.plot(t, data_homog.loc[t, 'median_filter'], marker='o', color='red', markersize=5)
for t in minimos.index:
plt.plot(t, data_homog.loc[t, 'median_filter'], marker='o', color='blue', markersize=5)
reconstr = data_homog['median_filter'].rolling(5, center=True).median()
for i, v in maximos.items():
reconstr.loc[i:] -= v
for i, v in minimos.items():
reconstr.loc[i:] -= v
reconstr.rolling(300, center=True).median().plot(ax=ax, color='magenta', lw=6, alpha=.4)
ax.xaxis.set_major_locator(mpd.HourLocator())
ax.xaxis.set_major_formatter(mpd.DateFormatter('%H:%M'))
In [87]:
data_homog['pendiente_5min'] = data_homog['median_filter'].rolling(300, center=True).apply(lambda x: x[-1] - x[0]).fillna(0)
pend_dia = data_homog.loc[data_homog.sun, 'pendiente_5min']
pend_noche = data_homog.loc[~data_homog.sun, 'pendiente_5min']
r_dia = (pend_dia.max() - pend_dia.min()) // 2
r_noche = (pend_noche.max() - pend_noche.min()) // 2
ax = pend_dia.hist(bins=r_dia, figsize=FS, log=True, alpha=.6, lw=0, color=tableau20[6])
pend_noche.hist(bins=r_noche, ax=ax, log=True, alpha=.7, lw=0, color=tableau20[8])
Out[87]:
In [91]:
sns.kdeplot(data_homog.loc[data_homog.sun, 'pendiente_5min'], data_homog.loc[data_homog.sun, 'median_filter'])
Out[91]:
In [127]:
#sns.kdeplot(data_homog.loc[~data_homog.sun, 'pendiente_5min'].values, data_homog.loc[~data_homog.sun, 'median_filter'].values)
ax = data_homog.loc[data_homog.sun, 'disp_r7'][data_homog['disp_r7'] > 10].hist(bins=100, figsize=FS, lw=0, alpha=.6, color='red')
data_homog.loc[~data_homog.sun, 'disp_r7'][data_homog['disp_r7'] > 10].hist(bins=100, lw=0, alpha=.6, ax=ax)
Out[127]:
In [ ]:
In [ ]:
In [170]:
dplot = data_homog.between_time('0:00', '23:59:59')
ax = dplot.drop(['poi', 'intervalo'], axis=1).iloc[::5].plot(figsize=FS)
ax.vlines(dplot[dplot['poi']].index, -600, 800, lw=1.5, color='y', alpha=.8)
(dplot['intervalo'] * 10).plot(ax=ax, lw=2, color='magenta', alpha=.8)
plt.grid('off')
ax.xaxis.set_major_locator(mpd.HourLocator())
ax.xaxis.set_major_formatter(mpd.DateFormatter('%H:%M'))
#data_homog.between_time('19:10', '20:00').plot(figsize=FS)
#data_homog.between_time('16:30', '17:30').plot(figsize=FS)
#data_homog['disp_r7'][data_homog['disp_r7'] > 20]
In [185]:
#print(dplot[dplot['poi']].head(7))
#dplot[(dplot.intervalo >= 0) & (dplot.intervalo <= 3)].plot(figsize=FS)
#dplot[dplot['poi']]
poi = dplot[dplot['poi']]
poi.set_index('intervalo').join(dplot.groupby('intervalo').median_filter.std().rename('std_mf')
).join(dplot.groupby('intervalo').median_filter.mean().rename('mean_mf'))
#dplot.groupby('intervalo').median_filter.mean()
Out[185]:
In [195]:
gb_int = dplot.tz_localize(None).reset_index().groupby('intervalo')
poi_data = gb_int.first().join(gb_int.median_filter.std().rename('std_mf')
).join(gb_int.median_filter.median().rename('median_mf')
).join(gb_int.median_filter.apply(lambda x: np.mean(x[:10])).rename('first_mf')
).join(gb_int.median_filter.apply(lambda x: np.mean(x[-10:])).rename('last_mf'))
poi_data
Out[195]:
In [194]:
plt.scatter(x=poi_data.median_mf, y=poi_data.median_filter.)
Out[194]:
In [60]:
#irrad_dif_dir
pysolar.solar.get_altitude_fast?
In [64]:
idx_sun = data_homog[data_homog.sun].index
irrad_dif_dir = pd.concat([pd.Series(idx_sun, index=idx_sun).iloc[::50].apply(lambda x: pysolar.util.diffuse_underclear(LAT, LONG, x)
).reindex(idx_sun).interpolate(),
pd.Series(idx_sun, index=idx_sun).iloc[::50].apply(lambda x: pysolar.radiation.get_radiation_direct(x, pysolar.solar.get_altitude_fast(LAT, LONG, x))
).reindex(idx_sun).interpolate()], axis=1)
irrad_dif_dir.loc[::20].between_time('8:00', '20:00').plot(figsize=FS)
Out[64]:
In [12]:
ldr = data_homog.ldr.values
@timeit('TOF', verbose=True)
def _TOF(arr):
n = 0
while n < 100:
total_variation = np.sum(arr[1:] - arr[:-1])
n += 1
return total_variation
@timeit('TOF_ne', verbose=True)
def _TOF_ne(arr):
n = 0
while n < 100:
#total_variation = ne.evaluate('sum()')
a, b = arr[1:], arr[:-1]
total_variation = ne.evaluate('sum(a - b)')
n += 1
return total_variation
_TOF(ldr)
_TOF_ne(ldr)
Out[12]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [3]:
data_solar.plot(figsize=(18, 8))
altura = data_solar.altitud[data_solar.altitud.abs() > 0]
azimut = data_solar.azimut[data_solar.altitud.abs() > 0]
In [ ]:
In [4]:
Image('http://www.monografias.com/trabajos82/energia-solar-fotovoltaica-y-sus-aplicaciones/image031.png')
Out[4]:
In [5]:
import math
def cosd(x):
return np.cos(np.radians(x))
def sind(x):
return np.sin(np.radians(x))
COS85_MIN = 0.0871557
def calc_cosenos(altitud_solar_deg, azimut_solar_deg, orientacion_N_deg, inclinacion_deg):
cenit_deg = 90 - altitud_solar_deg
azimut_efectivo = 180 + azimut_solar_deg - orientacion_N_deg
# Ángulos de incidencia
cos_theta = cosd(cenit_deg)*cosd(inclinacion_deg)+sind(cenit_deg)*sind(inclinacion_deg)*cosd(azimut_efectivo)
cos_incidencia = cos_theta.apply(lambda x: max(0, x))
cos_cenit_acotado = cosd(cenit_deg).apply(lambda x: max(COS85_MIN, x))
return cos_incidencia, cos_cenit_acotado
orientacion_N_deg, inclinacion_deg = 90 - 27.5, 90
out = pd.concat(calc_cosenos(altura, azimut, orientacion_N_deg, inclinacion_deg), axis=1)
ax = out.plot(figsize=(18,10))
orientacion_N_deg, inclinacion_deg = - 27.5, 90
out_2 = pd.concat(calc_cosenos(altura, azimut, orientacion_N_deg, inclinacion_deg), axis=1)
out_2.plot(ax=ax)
(altura / altura.max()).plot(ax=ax, lw=3)
Out[5]:
In [ ]:
In [ ]:
In [11]:
COS85_MIN = 0.0871557
def _cos_incidencia(data, orientacion_N_deg, inclinacion_deg=90,
alt='altitud', azi='azimut', irrad=None): # 'irradiacion_cs'):
cenit_deg = 90 - data[alt]
azimut_efectivo = 180 + data[azi] - orientacion_N_deg
# Ángulos de incidencia
cos_theta = cosd(cenit_deg)*cosd(inclinacion_deg)+sind(cenit_deg)*sind(inclinacion_deg)*cosd(azimut_efectivo)
cos_incidencia = cos_theta.apply(lambda x: max(0, x))
if irrad is not None:
if type(irrad) is str:
return cos_incidencia * data[irrad]
else:
return cos_incidencia * irrad
# cos_cenit_acotado = cosd(cenit_deg).apply(lambda x: max(COS85_MIN, x))
return cos_incidencia
In [17]:
data_solar_p = data_solar.copy()
labels = ['irradiacion_cs']
for ang_N in range(-180, 180, 18):
label = 'irrad_{}'.format(ang_N)
data_solar_p[label] = _cos_incidencia(data_solar_p, ang_N, irrad='irradiacion_cs', inclinacion_deg=90)
labels.append(label)
data_solar_p[labels].plot(figsize=(18, 12))
Out[17]:
In [94]:
from numpy import pi as PI
i_CONSTANTE_SOLAR_E_0 = 1367 # W/m2
i_COS_OBLICUIDAD_ELIPTICA_RAD = .917477 # cosd(23.44º)
i_SIN_OBLICUIDAD_ELIPTICA_RAD = .397789 # sind(23.44º)
i_DURACION_ANYO_SOLAR = 365.24 # W/m2
i_DIAS_DESDE_SOLST_INV_A_DIA_1 = 10.
i_DIAS_DESDE_DIA_1_A_PERIHELIO = 2.
i_SIGMA_STEFAN_BOLTZMANN = 5.67e-8 # W/(m^2·K^4)
i_VALOR_FORM_FACT_CENIT_PEREZ_K_RAD = 1.041
i_VALOR_COS_85_MINIMO = 0.0871557
i_DEG_A_RAD = PI / 180.
i_EMISIVIDAD_OPACOS_CTE = .9
i_RESISTENCIA_SUPERFICIAL_EXTERIOR = .04
i_PREC_COMPARA = .01
def get_julian_day_of_year(when):
jd = when.dayofyear
return jd + when.hour / 24. + when.minute / (24. * 60) + when.second / (24. * 3600)
def i_calcula_desviacion_eq_of_time_y_declinacion(dia_juliano):
ang_anual_unit_rad = 2. * PI / i_DURACION_ANYO_SOLAR
A = ang_anual_unit_rad * (dia_juliano + i_DIAS_DESDE_SOLST_INV_A_DIA_1)
B = A + 2. * .0167 * np.sin(ang_anual_unit_rad * (dia_juliano - i_DIAS_DESDE_DIA_1_A_PERIHELIO))
sin_declinacion = -i_SIN_OBLICUIDAD_ELIPTICA_RAD * np.cos(B)
cos_declinacion = np.cos(np.arcsin(sin_declinacion))
C = (A - np.arctan(np.tan(B) / i_COS_OBLICUIDAD_ELIPTICA_RAD)) / PI
EoT_opc = 720. * (C - np.floor(C + 0.5 + 0.0000001))
return sin_declinacion, cos_declinacion, EoT_opc
# Ecuation of time
def i_calcula_declinacion(dia_juliano):
ang_anual_unit_rad = 2. * PI / i_DURACION_ANYO_SOLAR
A = ang_anual_unit_rad * (dia_juliano + i_DIAS_DESDE_SOLST_INV_A_DIA_1)
B = A + 2. * .0167 * np.sin(ang_anual_unit_rad * (dia_juliano - i_DIAS_DESDE_DIA_1_A_PERIHELIO))
sin_declinacion = -i_SIN_OBLICUIDAD_ELIPTICA_RAD * np.cos(B)
cos_declinacion = np.cos(np.arcsin(sin_declinacion))
return sin_declinacion, cos_declinacion
def enhsol_calcula_altura_azimut_en_hora_solar(f_dia_anual, latitud_rad, con_correccion_eot):
def i_angulo_horario_de_f_dia_anual(f_dia):
hora_diaria = f_dia % 24. * 24.
return PI * (hora_diaria / 12. - 1.)
if con_correccion_eot:
sin_declinacion, cos_declinacion, correccion_EoT_min = i_calcula_desviacion_eq_of_time_y_declinacion(f_dia_anual)
else:
correccion_EoT_min = 0
sin_declinacion, cos_declinacion = i_calcula_declinacion(f_dia_anual)
angulo_horario = i_angulo_horario_de_f_dia_anual(f_dia_anual)
sin_alfa_s = np.sin(latitud_rad) * sin_declinacion + np.cos(latitud_rad) * cos_declinacion * np.cos(angulo_horario)
cos_gamma_s = (sin_declinacion - sin_alfa_s * np.sin(latitud_rad)) / (np.cos(latitud_rad) * np.sqrt(1 - sin_alfa_s**2))
altura = np.arcsin(sin_alfa_s) / i_DEG_A_RAD
if altura > 0.:
azimut = np.arccos(cos_gamma_s) / i_DEG_A_RAD
if (angulo_horario < 0.):
azimut = azimut - 180.
else:
azimut = 180. - azimut
if (con_correccion_eot == True):
azimut += correccion_EoT_min / 4.
else:
altura = 0.
azimut = 0.
return altura, azimut, 90. - altura # cenit_solar_deg
get_julian_day_of_year(ts)
print_red(i_calcula_desviacion_eq_of_time_y_declinacion(get_julian_day_of_year(ts)))
con_correccion_eot = True
enhsol_calcula_altura_azimut_en_hora_solar(get_julian_day_of_year(pd.Timestamp.now()), LAT * i_DEG_A_RAD, con_correccion_eot)
Out[94]:
In [34]:
ts = pd.Timestamp(day)
pysolar.solar.time.get_julian_solar_day(ts) % 365
Out[34]:
In [44]:
pysolar.solar.equation_of_time(get_julian_day_of_year(ts))
Out[44]:
In [60]:
compara_eot = pd.DataFrame([(i_calcula_desviacion_eq_of_time_y_declinacion(t)[2],
pysolar.solar.equation_of_time(t))
for t in range(365)], columns=['cype', 'pysolar'])
compara_eot.plot(figsize=(18, 10))
Out[60]:
In [100]:
jd_arr = data_solar.index.dayofyear + data_solar.index.hour / 24. + data_solar.index.minute / (24. * 60) + data_solar.index.second / (24. * 3600)
compara_alt_azi = pd.DataFrame([(*enhsol_calcula_altura_azimut_en_hora_solar(jd,
LAT * i_DEG_A_RAD, con_correccion_eot)[:2],
pysolar.solar.get_altitude(LAT, LONG, t),
pysolar.solar.get_azimuth(LAT, LONG, t)) for t, jd in zip(data_solar.index, jd_arr)],
columns=['alt_cype', 'azi_cype', 'alt', 'azi'], index=data_solar.index)
compara_alt_azi.plot(figsize=(18, 12))
Out[100]:
In [101]:
compara_alt_azi.apply(np.argmax)
Out[101]:
In [102]:
t = 240
i_calcula_desviacion_eq_of_time_y_declinacion(t)[2], pysolar.solar.equation_of_time(t)
Out[102]:
In [108]:
LONG / 360 * 24 * 60
Out[108]:
In [128]:
t
Out[128]:
In [119]:
# Prepara datos para comparación "fina":
#solar.iloc[::10].plot()
compara = solar[(solar.artif_level_max == 0) & (solar.altitud > 0)][['ldr', 'altitud', 'azimut', 'irradiacion_cs']].copy()
compara.plot()
Out[119]:
In [130]:
compara.altitud.argmax(), pysolar.solar.equation_of_time(pd.Timestamp('2016-08-27 14:05:06', tz=TZ).dayofyear)
Out[130]:
In [131]:
otro = get_solar_day(day=pd.Timestamp('2016-02-18'), step_minutes=1, tz=TZ)
otro.altitud.argmax(), pysolar.solar.equation_of_time(pd.Timestamp('2016-02-18 14:05:06', tz=TZ).dayofyear)
Out[131]:
In [133]:
long_minutes = LONG / 360 * (24 * 60)
long_minutes
Out[133]:
In [135]:
(pysolar.solar.equation_of_time(pd.Timestamp('2016-02-18 14:05:06', tz=TZ).dayofyear),
pysolar.solar.equation_of_time(get_julian_day_of_year(pd.Timestamp('2016-02-18 14:05:06', tz=TZ))),
pysolar.solar.equation_of_time(pd.Timestamp('2016-02-19 14:05:06', tz=TZ).dayofyear))
Out[135]:
In [246]:
data_solar_p = compara.iloc[::100].copy()
labels = []
for ang_N in range(-180, 180, 5):
label = 'irrad_{}'.format(ang_N)
data_solar_p[label] = _cos_incidencia(data_solar_p, ang_N, irrad='irradiacion_cs', inclinacion_deg=90)
labels.append(label)
ax = data_solar_p[['ldr', 'irradiacion_cs']].plot(figsize=(18, 12), lw=3)
for c in labels:
data_solar_p[c].plot(ax=ax, lw=.75, legend=c)
In [150]:
data_solar_p = compara.iloc[::100].copy()
labels = []
for ang_N in range(-180, 180, 30):
label = 'irrad_{}'.format(ang_N)
data_solar_p[label] = _cos_incidencia(data_solar_p, ang_N, irrad='irradiacion_cs', inclinacion_deg=10)
labels.append(label)
ax = data_solar_p[['ldr', 'irradiacion_cs']].plot(figsize=(18, 12), lw=3)
for c in labels:
data_solar_p[c].plot(ax=ax, lw=.75, legend=c)
In [214]:
data_solar_p = compara.iloc[::100].copy()
labels = []
for ang_N in range(-180, 180, 5):
label = 'irrad_{}'.format(ang_N)
data_solar_p[label] = _cos_incidencia(data_solar_p, ang_N, irrad='irradiacion_cs', inclinacion_deg=45)
labels.append(label)
ax = data_solar_p[['ldr', 'irradiacion_cs']].plot(figsize=(18, 12), lw=3)
for c in labels:
data_solar_p[c].plot(ax=ax, lw=.75, legend=c)
In [247]:
import statsmodels.api as sm
X, Y = data_solar_p['ldr'], data_solar_p.drop(['altitud', 'azimut', 'ldr'], axis=1)
model = sm.OLS(X, Y).fit()
model.summary()
(model.params * Y).sum(axis=1).plot(figsize=(18, 10), color='darkred')
X.plot(lw=3, alpha=.7, color='darkgreen')
model.pvalues[model.pvalues > .5]
Out[247]:
In [248]:
(model.params[model.params > 0] * Y[model.params[model.params > 0].index]).sum(axis=1).plot()
Out[248]:
In [278]:
from enerpi.base import timeit
@timeit('separa_ldr_artificial_natural', verbose=True)
def _separa_ldr_artificial_natural(data_day, resample_inicial='2s', delta_roll_threshold=100):
def _f_roll(x):
delta_pos = np.sum(x[x > 0])
delta_neg = np.sum(x[x < 0])
return delta_pos + delta_neg
def _last_on(x):
if (np.sum(x) == 1.) & (x[-1] == 1.):
return True
return False
# TODO Resolver NaN's
data_analysis = pd.DataFrame(data_day['ldr'].resample(resample_inicial).mean().interpolate().dropna())
data_analysis['delta_roll'] = data_analysis.ldr.diff().fillna(0).rolling(3, center=True).apply(_f_roll).fillna(0)
data_analysis['ch_state'] = 0
data_analysis.loc[data_analysis[(data_analysis['delta_roll'] > delta_roll_threshold
).rolling(3).apply(_last_on).fillna(False) > 0].index, 'ch_state'] = 1
data_analysis.loc[data_analysis[(data_analysis['delta_roll'] < -delta_roll_threshold
).rolling(3).apply(_last_on).fillna(False) > 0].index, 'ch_state'] = -1
cambios = data_analysis[data_analysis['ch_state'] != 0]
print_red(cambios)
data_analysis['artif_level_max'] = 0
apagados = data_analysis[data_analysis['ch_state'] == -1].index
for i_enc in data_analysis[data_analysis['ch_state'] == 1].index:
apagado = apagados[apagados > i_enc]
if len(apagado) > 0:
apagado = apagado[0]
data_analysis.loc[i_enc:apagado, 'artif_level_max'] = data_analysis.ldr.loc[i_enc:apagado].max()
print_red('Encendido L={:.0f}, SEGS={:.0f}'.format(data_analysis.ldr.loc[i_enc:apagado].max(),
(apagado - i_enc).total_seconds()))
# Reconstrucción LDR natural:
hay_artificial = data_analysis.artif_level_max > 0
index_art = data_analysis[hay_artificial.shift(-2) | hay_artificial.shift(2)].index
rs_natural = data_analysis.ldr.drop(index_art).resample('2min')
data_analysis_simple = pd.concat([rs_natural.max().rename('nat_max').interpolate(),
rs_natural.min().rename('nat_min').interpolate(),
rs_natural.mean().rename('nat_mean').interpolate()], axis=1)
# data_analysis_simple = data_analysis_simple.join(data_analysis[['altitud', 'azimut']])
return data_analysis, data_analysis_simple
day = '2016-09-09'
d_sample = LDR.loc[day:day]
d_sample, d_sample_simple = separa_ldr_artificial_natural(d_sample,
resample_inicial='2s', delta_roll_threshold=100)
d_sample.iloc[::100].plot()
Out[278]:
In [279]:
d_sample.iloc[::100].delta_roll.plot()
Out[279]:
In [304]:
def _roll_band(x):
ab = np.max(x) - np.min(x)
#if ab > 15:
# return ab
#return 0
return ab / np.mean(x)
d_sample.ldr.rolling(15, center=True).apply(_roll_band).plot()
Out[304]:
In [305]:
#((d_sample.ldr).rolling(15).mean().pct_change()).plot()
(d_sample.ldr.rolling(15, center=True).apply(_roll_band) * d_sample.ldr).plot()
Out[305]:
In [297]:
((d_sample.ldr).rolling(15).mean().diff()).plot()
Out[297]:
In [222]:
cols_pos = model.params[model.params > 0].index
model_p = sm.OLS(X, Y[cols_pos]).fit()
model_p.summary()
Out[222]:
In [229]:
model_p.pvalues[model_p.pvalues > .5]
Out[229]:
In [221]:
(model_p.params * Y[cols_pos]).sum(axis=1).plot()
X.plot()
Out[221]:
In [172]:
dir(scipy.optimize)
scipy.optimize.curve_fit?
In [220]:
sol_v, rnorm = scipy.optimize.nnls(Y, X)
ax = (sol_v * Y).sum(axis=1).plot()
(sol_v * Y).plot(ax=ax)
X.plot()
Out[220]:
In [197]:
def local_maxima_minima(serie, window_delta='2h', verbose=True):
def _is_local_min(x):
if np.argmin(x) == len(x) // 2:
return True
return False
def _is_local_max(x):
if np.argmax(x) == len(x) // 2:
return True
return False
roll_w = int(pd.Timedelta(window_delta) / serie.index.freq)
if roll_w % 2 == 0:
roll_w += 1
if verbose:
print_secc('LOCAL MAX-MIN de {} valores. Centered window of {} samples (for ∆T={}, freq={})'
.format(len(X), roll_w, window_delta, serie.index.freq))
maximos = serie.rolling(roll_w, center=True).apply(_is_local_max).fillna(0).astype(bool)
minimos = serie.rolling(roll_w, center=True).apply(_is_local_min).fillna(0).astype(bool)
if verbose:
print_info(serie[maximos])
print_cyan(serie[minimos])
return maximos, minimos
maximos, minimos = local_maxima_minima(X, window_delta='2h', verbose=True)
ax = X.plot(figsize=(18,10))
X[maximos].plot(ax=ax, lw=0, marker='o', color='red')
X[minimos].plot(ax=ax, lw=0, marker='o', color='blue')
#ax.set_xlim((X[maximos].index[0], X[maximos].index[-1]))
Out[197]:
In [ ]:
In [ ]:
In [ ]: